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Abstract 

A two-dimensional unstructured Navier-Stokes code is utilized for computing the flow around multielement airfoil configurations. 
Comparisons are shown for a landing configuration with an advanced-technology flap. Grid convergence studies are conducted to assess 
inaccuracies caused by inadequate grid resolution. Although adequate resolution is obtained for determining the pressure distributions, 
further refinement is needed to sufficiently resolve the velocity profiles at high angles of attack. 

For the advanced flap configuration, comparisons of pressure distributions and lift are made with experimental data. Here, two flap 
riggings and two Reynolds numbers are considered. In general, the trends caused by variations in these quantities are well predicted by 
the computations, although the angle of attack for maximum lift is overpredicted. 


Introduction 

The goal of a high-lift system is to generate as much lift 
as possible without separating the flow [1], Without external 
devices such as wall suction, the most effective way to achieve 
this goal is through the use of multiple elements to manipulate the 
inviscid pressure distribution to reduce the pressure rise over each 
element [1], [2], However, the presence of multiple elements 
seriously complicates analysis procedures because of important 
and often complex interactions between the individual elements. 
While inviscid analysis can be accomplished in minutes with panel 
methods or unstructured-grid Euler solvers, it is necessary to use 
viscous techniques to accurately predict the flows about these 
configurations. The reason for this is that although the tailoring 
of the flow field to prevent separation is largely achieved through 
circulation interactions between the elements, many viscous effects 
can have large influences on the pressure distributions. While 
these include obvious effects such as displacement thickness and 
separation on the surfaces, wake interactions between forward and 
aft elements as well as flow reversal off the surface can contribute 
significantly in determining the overall performance of the high- 
lift system [3], [4], [5], 

For computations on multielement airfoils, unstructured grid 
methods may offer a good alternative to more traditional methods 
of analysis. This is due in part to the decreased time required to 
generate grids over complicated geometries. Also, unstructured 
grids offer the potential to adapt the grid to improve the accuracy 
of the computation without incurring the penalties associated with 
global refinement. However, despite the advantages of unstruc- 
tured grids, they are typically much slower than structured grid 
solvers. Also, the ability to obtain solutions through local adapta- 
tion that are comparable to those obtained through global refine- 
ment remains an area where further work is required [ 6 ]. Although 
work remains to fully realize their potential, much progress has 
been reported in computing viscous flows on unstructured grids 
(See for example [7], [ 8 ], [9], [10]). 

The purpose of this study is to present computational results 
obtained with a particular unstructured grid method [ 11 ], [ 12 ]that 
has been applied to several flows over multielement airfoils. Com- 
parisons between computational results and experimental data are 
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made to assess the effectiveness of the present code, to aid in de- 
termining future directions, and to provide useful comparisons for 
other researchers working in this field. 

Symbols 


a angle of attack 

b span of wind tunnel model 

c reference chord taken to be the chord of the 

undeflected airfoil 

C.'i lift coefficient 

C p pressure coefficient 

CFL Courant-Friedichs-Lewy number 

.1/ free stream Mach number 

q magnitude of velocity 

7 magnitude of velocity in freestream 

Re Reynolds number 

u velocity component in direction of 

surface-tangent vector 

x, y, z Cartesian coordinates 

y+ turbulent boundary-layer parameter 

a angle of attack 

</ chordwise location on airfoil (referenced to 

undeflected position) 

Computational Method 


The computational method used in this study is a node-based, 
implicit, unstructured, upwind flow solver described in reference 
[1 1], In this code, the discretization of the convective and viscous 
terms is handled similarly to the method of reference [ 8 ], The 
inviscid fluxes are obtained using Roe’s approximate Riemann 
solver [13]; the viscous terms are evaluated with a Galerkin-type 
approximation that results in a central-difference formulation for 
these terms. Two different turbulence models are presently utilized 
in the code. These include both the Baldwin-Barth [14] model and 
the Spalart-Allmaras [15] model. At each time step, the equation 
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for the turbulent viscosity is solved separately from the flow 
equations, which results in a loosely coupled solution process that 
allows for a relatively easy interchange of other turbulence models. 
Although both turbulence models have been used extensively with 
good success, the present study reports only results obtained with 
the Spalart-Allmaras model. 

Experimental Data 


All experimental data used in the present work have been 
obtained in the Low Turbulence Pressure Tunnel (LTPT) located 
at the NASA Langley Research Center [16]. The tunnel is a 
single return, closed-throat wind tunnel that obtains high Reynolds 
numbers by operating at pressures up to 10 atm. The test section 
is 3 ft wide by 7.5 ft high by 7.5 ft long. Side wall boundary 
suction is applied to promote two-dimensional flow [17]. 

Lift and moment measurements are obtained by using both 
a force balance and an integration of surface pressures; drag is 
obtained from a wake survey with a five-hole probe. The accuracy 
of the lift force is approximately ±0.5 percent when obtained 
from the balance [16]. and the lift coefficient is estimated to be 
within ±0.03 when obtained from pressure integration. The drag 
coefficient is estimated to be accurate to within ±0.001 [18] for 
attached flows. Pressure coefficient distributions are obtained from 
pressure orifices located along the model and are accurate to within 
about ±0.030. 

For the calculations that follow, comparisons will be made 
with experiments obtained from two different data sets and the 
presentation of results is organized accordingly. A brief descrip- 
tion is given below for each data set. 

The test data is the result of a cooperative experimental pro- 
gram between the Douglas Aircraft Company and the NASA Lan- 
gley Research Center and is reported in references [19] and [20], 
Test Reynolds numbers varied between 5, 9, and 16 million. The 
angle of attack included a range of approximately —4° through 
23°. The tests have been conducted without forced boundary- 
layer transition. The overall geometry, which is shown in figure 
1, is a three-element configuration based on an 11.55 percent thick 
supercritical airfoil. The slat and flap chord ratios are 14.48 per- 
cent and 30 percent, respectively, based on the airfoil chord for 
the undeflected position. 

For the current study, the deflections of both the slat and the 
flap are set at 30°. and two different flap riggings are considered. 
A “rigging” refers to a combination of gap and overhang settings 
as defined in figure 2, and a specific rigging is assigned a letter 
designation. For the first configuration, denoted as 30P-30N (slat 
deflection 30°, slat rigging P, flap deflection 30°, flap rigging N), 
the flap overhang is 0.25 percent of the undeflected airfoil and the 
flap gap is 1.7 percent. For the second configuration, which is 
denoted as 30P-30AG, both the gap and overhang are 1 percent. 
In figure 3, a more detailed view of the flap riggings for the two 
configurations is shown. 

The data for these configurations have been obtained from 
two separate tunnel entries. Because the first test considered both 
three- and four-element airfoils, the flap was constructed in two 
pieces, which were then assembled on site. During these tests, 
force and moment data were obtained for many flap riggings, 
including those discussed here. From these experiments, it was 
found that the 30P-30N configuration exhibited a slightly higher 
( 'i liirl , than did the 30P-30AG. A single-segment flap was then 
constructed and the 30P-30N geometry was studied in more detail 
in a subsequent test. During this test, more detailed data, such 


as velocity profiles, were obtained and are presented in reference 

[20], 

In order to make meaningful interpretations between compu- 
tational and experiment results, an indication of the two dimen- 
sionality of the flowfield is necessary. Figure 4 shows experimen- 
tal pressure distributions at several angles of attack obtained at two 
different spanwise locations on the model. The z/b = 0.5 location 
corresponds to the centerline of the tunnel which has 146 pressure 
taps. In order to better assess the two dimensionality of the flow, 
several pressure taps are also present at z/b = 0.77, which is ap- 
proximately midway between the centerline of the model and the 
wind tunnel wall. As seen in the figure, excellent two-dimensional 
flow is maintained at an angle of attack of 16.2° (uncorrected). 
Slight three-dimensional effects are present at 21.31° and three 
dimensionality is clearly indicated at 22.25°. 

Results 

The results of the comparisons are presented below. For the 
computations, the CFL number has been ramped linearly from 20 
to 100 over 100 iterations. The lift coefficient is constant to within 
about 0.1 percent of the total lift (fourth significant digit) and the 
pressure distributions and velocity profiles show no differences 
over 500 iterations when plotted. 

All grids used in the present study have been generated 
using the grid generation procedure described in reference [21]; a 
sample is shown in figure 5. In this procedure, structured grids 
are first generated around individual components. These grids 
are then used to define a cloud of points that is triangulated 
with a stretched-Delaunay triangulation procedure to establish 
the connectivity relationships. Although the aspect ratios of the 
triangles near the surface are generally very large because of the 
extremely small spacing required at the wall, grids generated in 
this manner tend to have relatively few cells with large angles 
that can negatively compromise accuracy. For the grids used here, 
fewer than 2 percent of the angles are greater than 120°. 

Comparisons between computed results and experimental 
data are presented below for the three-element advanced flap 
configurations (30P-30N and 30P-30AG) discussed above. The 
first set of results is for the 30P-30N configuration at several angles 
of attack with a Mach number of 0.2 and a Reynolds number of 9 
million. Simultaneously presented are the results of varying grid 
densities on computed pressure distributions and velocity profiles, 
which are used to ascertain the level of numerical errors in the 
computations. 

For these studies, three grids have been utilized. The first 
grid, which will be referred to as the coarse grid, consists of 
22,491 nodes, 524 of which lie on the surfaces of the elements. 
The spacing at the wall for this mesh is 4 x 10 ~ 6 normalized to 
the chord length of the airfoil in the undeflected position. This 
spacing yields a of less than 2 for the point next to the wall 
at the rear of a unit-length flat plate. Finer grids are obtained by 
simultaneously increasing the number of points in each direction 
to obtain as close to a uniform refinement as possible. In the 
current study, the number of points in each direction is increased 
with each refinement by a factor of roughly \/2 so that the total 
number of nodes with each refinement is approximately doubled 
over the previous mesh. 

With this procedure, the second mesh in this family of grids 
contains 49,596 nodes with 806 points on the surfaces. This mesh 
will be referred to in future discussions as the baseline mesh. For 
this mesh, the spacing at the wall is about 2 x 10 “ 6 and results in 
a of less than 1, according to flat plate estimates. In a similar 
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fashion, another finer grid is generated that has 87,783 nodes with 
1036 points on the solid surfaces; it is referred to here as the fine 
grid. 

Pressure distributions on all grids are shown in figures 6 
through 9 for a = —0.03°, a = 8.23°, a = 16.30°, and 
a = 22.36°, respectively. As seen, the variation in the pressure 
distribution between the coarse grid and the other two grids is 
relatively small up to a = 22.36°. However, at this angle of 
attack, the loading on the flap is noticeably less for the coarse 
grid than for the baseline and fine grids, which continue to yield 
very similar results. Also, an obvious discrepancy exists in the 
pressure distribution on the slat for an angle of attack of 8.23°. 
The cause of this is unknown, but could be attributable to a 
number of sources including wind tunnel wall corrections or an 
inaccurate representation of the separated flow under the cove 
which leads to higher circulation around this element. Although 
not shown, it is interesting to note that excellent agreement is 
obtained on the slat, as well as the other elements, when the 
point vortex correction in the far field for the computations is 
removed. However, numerical studies for the Euler equations with 
this same three-element airfoil indicate that the presence of the 
point vortex significantly decreases the dependence of the solution 
on the placement of the far-field boundary as expected. Therefore, 
the improved agreement obtained without the vortex is evidently 
fortuitous. Note also that the pressure distributions on both the 
main element and the flap indicate slightly higher lift than the 
experiment, which could also account for a higher circulation on 
the slat. 

Computed velocity profiles at several locations along the 
main element and the flap are shown in figure 10 for an angle 
of attack of 22.36°. Also shown is an illustration that indicates 
the locations on each element where the data are obtained. Note 
that these locations are referenced to the airfoil coordinates in 
the undeflected position. Results are shown for this angle of 
attack because they indicate the most variation between the coarse, 
baseline, and fine grids. As seen, the agreement between the 
baseline grid and the fine grid is reasonably good on the main 
element but the coarse grid is clearly inaccurate. Furthermore, 
the trend with grid refinement is to decrease the boundary-layer 
thickness on the main element. The increased thickness of the 
boundary layer with the coarse grid is likely to be responsible 
for the decreased flap loading obtained on this grid and seen in 
figure 9. 

On the flap, the major difference between the baseline grid 
results and the fine grid results is in the enhanced resolution of 
the slat wake on the finer grid. This wake is very apparent at the 
station immediately downstream of the main element (// = 0.72) 
but quickly dissipates so that its presence is barely detectable 
at the location towards the back of the flap (jj = 0.92). The 
difference in profiles between the baseline and fine grids indicates 
that further refinement is necessary to accurately resolve these 
details. Although not shown, as the angle of attack is reduced, 
the difference in the profiles decreases so that the baseline grid 
and fine grid give essentially identical results at an angle of attack 
of -0.03°. 

An additional study of grid effects has been conducted and 
is presented in figure 1 1 . For this study, a grid very similar to the 
baseline grid has been generated with the spacing of the grid points 
next to the wall based on obtaining a y + « 10 instead of // + ss 1. 
Note that the spacing is determined based on estimates from a 
flat plate at the Reynolds number in question (Re = 9 x 10 6 ). 
In order to obtain y + « 1, the required spacing at the wall is 
approximately 2 x 10 “ 6 ; however, y + « 10 allows spacing an 


order of magnitude larger. Although not shown here, the y + 
obtained from actual computations on the baseline grid at an angle 
of attack of 22.36° is approximately one over the first 20 percent 
of the airfoil and then drops to about 1/2 afterwards. Values of y + 
for the second mesh are slightly above 10 for the first 20 percent 
and drop to around 5 for most of the remainder of the element. 
Because the sublayer for a turbulent boundary layer extends to a 
y + of approximately 10, essentially no points exist in this region 
over the first 20 percent of the airfoil for this mesh. 

Figure 11 shows that inadequate spacing at the wall drasti- 
cally affects the pressure distribution on the flap; hence, the other 
elements are effected as well. Inspection of the velocity profiles 
(not shown) reveals that inadequate spacing near the wall leads 
to an artificial thickening of the computed boundary layer on the 
main element. As discussed in reference [3], a thick boundary 
layer on the main element acts to suppress the loading on the flap. 
Therefore, the effect of the artificially thickened boundary layer is 
to artificially suppress the loading on the flap. 

Numerical experiments indicate that inadequate wall spacing 
has little effect on the pressure distributions at lower angles of 
attack, such as 8.23°. It is primarily at higher angles where the 
wall spacing has been observed to be critical. Note that for single 
element airfoils, the effect of wall spacing is not as dramatic as 
it is for multielement configurations because the wake does not 
impinge on aft elements. Also, achieving a y + of 1 is not a 
necessary requirement; typically, values similar to those used for 
algebraic turbulence models (</ + « 3) should be sufficient [14], 

[15]. 

A summary of computed and experimental lift coefficients 
at Reynolds numbers of both 5 and 9 million are shown in figure 
12. Here, the lift versus angle of attack is shown for the full 
configuration, as well as for the individual elements. The lift 
for each of the individual elements is obtained from pressure 
integration and is not corrected for wind tunnel wall effects. 
The total lift, on the other hand, is also computed from pressure 
integration, but has been corrected for wall interference. 

As seen in the figure, the lift agreement is good for the un- 
corrected data on each element whereas the total lift (corrected) for 
the configuration is overpredicted in the computations. Although 
not shown, improved agreement between the computations and the 
experiment is obtained by using the lift from the force balance be- 
cause it is slightly higher. For the computations, the angle of attack 
for maximum lift is not accurately predicted for either Reynolds 
number; however, the overall trend is well represented. The ex- 
periment and computations both obtain higher lift for a Reynolds 
number of 9 million over that of 5 million, and the main element 
begins to lose lift before either of the other elements. 

A comparison of velocity profiles for Reynolds numbers of 
5- and 9 million are shown in figure 1 3 for a = 16.3°. Only small 
differences are seen over the main element, although the boundary 
layer at 9 million appears to be slightly thinner than at 5 million. 
Over the flap, more differences are apparent; the computations 
at the lower Reynolds number show lower velocities than at a 
Reynolds number of 9 million. Although the results from the grid 
convergence studies indicate that more refinement is required to 
adequately resolve the wake emanating from the slat which persists 
to the back of the flap, the overall trends in the velocity profiles 
with variations in Reynolds number are well captured. 

The last case considered from this data set is a comparison 
of the computations and experiment between the 30P-30N and 
the 30P-30AG configurations. Recall from the description of the 
wind tunnel tests that these configurations differ only in the flap 
rigging, as shown in figure 3. Because the tests for the 30P-30AG 
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configuration were conducted with only the two-segment flap, the 
data shown below for both the 30P-30N and the 30P-30AG are 
from this same test. 

A comparison of the lift coefficients versus angle of attack for 
these configurations is shown in figure 14 for a Reynolds number 
of 9 million. In both the computations and the experiment, the 
30P-30N configuration attains higher lift coefficients than the 30P- 
30AG. Note that the agreement between the computed lift and 
the experimental data for the 30P-30N configuration is somewhat 
improved over that shown previously in figure 12. This is because 
the data for the two-segment flap is used for the current figure. 
Careful inspection of figures 12 and 14 indicates that differences in 
the experimentally obtained lift between the one-and two-segment 
flap appear to start at about 20°. This difference may be due to 
slight differences in the model geometry or because the flow may 
no longer be two dimensional at this angle of attack. 

A comparison of computed and experimental pressure distri- 
butions is shown in figure 15 for an angle of attack of 16.3°. The 
agreement between the computations and experiment is good for 
all elements. Although little difference due to the flap rigging is 
apparent on the slat and the main element, the flap shows a higher 
suction peak for the 30P-30N configuration than for 30P-30AG. 

Conclusions and Recommendations 

A two-dimensional unstructured Navier-Stokes code has 
been utilized for computing the flow around multielement airfoil 
configurations. Comparisons are shown for a landing configura- 
tion with an advanced-technology flap for angles of attack up to 
the maximum-lift condition. A systematic grid convergence study 
has been conducted to assess the inaccuracies in the computations 
caused by inadequate grid resolution. Below maximum lift, pres- 
sure distributions are adequately resolved by using approximately 
50,000 nodes. However, at high angles of attack, further grid re- 
finement is required to obtain suitable levels of grid convergence 
for velocity profiles. This could be achieved by continuing to re- 
fine the mesh in a systematic manner, or possibly through the use 
of adaptive gridding or higher-order methods. The grid studies 
further indicate that care must be taken in obtaining accurate reso- 
lution of the wall boundary layers on upstream elements by using 
sufficiently small spacing of grid points. The use of a grid with 
inadequate wall spacing (y + = 0(10)) results in an artificially 
thick boundary layer on the main element that severely effects the 
loading on the flap and hence the entire configuration. 

Comparisons for the advanced flap configuration between 
computed and experimental pressure distributions are made for 
two flap riggings. In addition, lift coefficients and velocity profiles 
are compared for Reynolds numbers of 5- and 9 million. In 
general, the trends due to variations in rigging and Reynolds 
numbers are predicted well by the computations although the angle 
of attack for maximum lift is overpredicted. 
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Figures 



Figure 1. Geometry for three-element airfoil 30P-30AG. 



x, in 



Figure 2. Definition of gap and overhang for flap. 

a=16.20° 


Figure 3. Differences in flap rigging for 
the 30P-30N and 30P-30AG configurations. 

a=21.31° a=22.25° 
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Figure 4. Experimental pressure distributions at two spanwise locations for several angles of attack. 



Figure 5. View of sample unstructured grid for three-element airfoil. 


5 


ANDERSON ET AL.: NAVIER-STOKES COMPUTATIONS 


Slat Main Flap 





Figure 6. Computed pressure distributions on the coarse, baseline, and fine grids for 30P-30N with Moo =0.2, R e =9xl0 6 ,o = —0.03° 



Figure 7. Computed pressure distributions on the coarse, baseline, and fine grids for 30P-30N with M = 0.2, R e = 9 x 10 6 , a = 8.23° 
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Figure 8. Computed pressure distributions on the coarse, baseline, and fine grids for 30P-30N with .1/ = 0.2, R e = 9 x 10 6 , a = 16.30° 



Figure 9. Computed pressure distributions on the coarse, baseline, and fine grids for 30P-30N with Moo = 0.2, R e = 9xl0 6 ,o = 22.36° 
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Figure 14. Comparison of computational and experimental lift for the 
30P-30N and 30P-30AG configurations with \l = 0.2. R e = 9 x 10 6 . 
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Figure 15. Comparison of computational and experimental pressure distributions for the 
30P-30N and 30P-30AG configurations with M, x . = 0.2, R e = 9 x 10 6 , a = 16.3°. 
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